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ABSTRACT 


The  absence  of  body  force  and  the  natural  orientations  of  two  phase- 
gas-liquid  systems  in  space  produce  major  differences  between  liquid  fuel 
transfer  in  space  and  that  on  earth.  A vortex-induced,  liquid  handling 
process  adapted  to  two-phase  fluids  in  zero  gravity  has  been  analyzed  for  a 
selected  range  of  fluid  and  flow  parameters.  This  transfer  process  has  been 
divided  in  two  stages:  the  initial  "spin-up”  stage  and  the  liquid  "pump- 
out" stage  from  the  "sender"  tank.  The  initial  spin-up  stage  is  established 
by  tangential  fluid  injection.  The  model  is  based  on  a two-phase,  two-fluid 
continuum  with  several  phase  interactions  - namely  fluid  drag  and  pertinent 
virtual  mass  effects. 

A computer  program  was  developed  to  study  the  fluid  dynamical  behavior 
of  two-phase  fluids  in  the  sender  tank.  Several  examples  are  given  to 
demonstrate  the  products  of  the  program.  The  interesting  and  plausible 
results  indicate  the  simplified  two-phase,  two-fluid  model  is  a very  useful 
one.  However,  detailed  evaluations  of  the  results  are  only  possible  with 
accurate  local  measurements . These  verifications  should  be  devised  and 
carried  out  to  check  both  the  formulations  selected  for  phase  interactions 
as  well  as  the  results  predicted. 

Key  Words:  Fuel  transfer  in  zero  gravity,  gas-liquid  separation,  numerical 
modeling,  rotational  flows,  two-phase  flows,  vortex  motions, 


vortex-induced  phase  separations. 
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1 . Introduction 


This  work  is  mainly  concerned  with  liquid  fuel  transfer  in  space,  i.e., 
in  zero  gravity.  A number  of  concepts  for  liquid  transfer  between  "sender" 
and  "receiver"  tanks  in  a zero-gravity  environment  have  been  proposed  [1], 
All  of  these  are  concerned  in  part  with  the  problem  of  collecting  the  liquid 
phase  and  feeding  it  to  a pumping  unit  while  providing  for  the  return  of 
gases  from  the  receiver  tank  to  the  sender  tank.  The  major  difference 
between  the  fuel  transfer  in  space  and  that  on  earth  is  the  absence  of  body 
forces  in  space.  Thus,  there  is  no  natural  orientation  of  the  liquid  and 
vapor  in  a tank. 

The  transfer  of  fuel  by  pumping  requires  separation  of  the  liquid  and 
gas.  The  basic  requirement  is  to  orient/collect  the  liquid  to  be 
transferred  in  the  sender  tank  and  to  provide  the  vent  system  in  the 
receiving  tank  to  prevent  excessive  liquid  loss.  Depending  on  different 
situations,  there  are  several  techniques  used  to  separate  liquid  from  gas  in 
order  to  optimize  the  transfer  process. 

Positive  displacement,  linear  acceleration,  tumble  acceleration,  roll 
acceleration,  centrifugal  fluid  motion,  dielectrophoresis , acoustic,  and 
magnetic  forces,  and  surface  tension,  etc.  are  the  most  common  techniques 
referred  [1],  Each  technique  has  its  own  shortcoming.  For  example, 
rotation  of  the  entire  vehicle  and  receiver  was  not  considered  practical  due 
to  adverse  dynamic  effects  and  changing  center  of  gravity  while 
transferring.  Rotation  of  the  tankage  within  the  vehicle  was  not  desirable 
due  to  the  requirement  for  stationary  to  rotational  connections.  An 
alternative  to  the  roll  acceleration  method  is  to  create  a centrifugal 
motion  in  the  fluid  while  holding  the  vehicle  system  stable. 
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Fluid  rotation  can  be  provided  by  various  mechanical  systems  or  by 
fluid  injection.  The  main  disadvantages  of  mechanical  systems  are  their 
requirements  for  motor  drives  and  potentially  higher  residual  weights  (e.g., 
long  diameter  paddles  that  would  sweep  liquid  away  from  the  wall  while  short 
diameter  paddles  might  not  be  effective  for  partially  filled  tanks). 
Therefore,  in  this  report,  fluid  rotation  created  by  fluid  injection  is 
considered . 

The  centrifugal  forces  created  by  vortex  motion  have  been  used  for 
various  separation  applications.  The  hydrocyclone  used  to  separate  solid 
particles  via  specific  gravity  differences  is  a very  economical  hydraulic 
device  with  no  moving  parts.  However,  some  basic  questions  are: 

(1)  Can  a vortex  motion  be  used  practically  to  orientate  the  fuel 
location  for  fuel  transfer  in  zero  gravity? 

(2)  If  so,  can  it  do  the  job  more  effectively  than  other  methods? 

To  answer  these  questions,  one  requires  a more  complete  understanding 
of  the  hydraulic  character istics  of  two-phase  vortex  motion.  Normally  the 
residuals  are  dependent  on  the  volumetric,  pump-out  flowrate  from  the 
tanker.  Intuition  suggests  that  one  should  set  a desirable  acceleration 
level  for  pump-out,  select  the  optimum  flow  rate  and  control  this  flow  rate 
toward  the  end  of  the  transfer  time.  The  basic  approach  is  to  determine  at 
what  rate  the  fluid  must  be  rotated  to  insure  liquid  at  the  wall.  Thus,  it 
is  desirable  to  properly  describe  the  pertinent  flow  processes  and  to 
produce  the  proper  set  of  equations  that  can  be  manipulated  and  solved  to 
give  the  necessary  answers  for  design  purposes. 

To  help  in  answering  the  above  questions  and  to  provide  design 
guidance,  we  have  developed  a physical/mathematical  model  for  the  fluid 
dynamical  behavior  of  liquid-gas  flow  in  a cylindrical  tank  at  zero-gravity. 


The  general  objectives  of  the  study  are:  1)  to  develop  a computer  model  for 
a two-phase  vortex  motion  at  zero-gravity;  2)  to  study  fluid  dynamical 
behavior  of  the  two  phase  in  the  sender  tank;  3)  to  provide  a mathematical 
framework  containing  the  essential  physical  mechanisms  on  vortex  induced 
phase-separation;  4)  to  provide  quick/easy  sensitivity  tests  on  various 
parameters;  and  5)  to  provide  guidance  for  the  design  and  testing  of  MIST 
(Mechanically  Induced  Settling  Technology)  transfer  systems,  and  their 
scale-up  to  appropriate  sizes. 

2.  Transfer  Concept  and  Mode  of  Operation 

The  MIST  Transfer  Concept  recently  proposed  by  KSC-NASA  uses  fluid 
induced  rotation  (vortex  motion)  to  establish  an  annular  region  of  the 
liquid  phase  adjacent  to  the  walls  of  a cylindrical  tank. 

Figure  1 presents  a simplified  diagram  of  the  MIST  Transfer  System. 
Initially,  the  gas  and  liquid  phases  in  the  sender  tank  are  arbitrarily 
distributed.  During  the  start-up  stage,  the  liquid  will  be  centrifuged 
outwardly  while  the  gas  will  move  radially  inward.  The  time  required  to 
reach  continuous  liquid  transfer  will  vary  depending  upon  the  initial 
distribution  of  phases.  To  provide  for  optimum  operation  of  the  MIST 
transfer  system  and  reliable  design  and  scale-up,  the  fluid  dynamical 
behavior  of  the  fluids  in  the  sender  tank  must  be  properly  modeled.  In  this 
work,  the  objective  was  to  estimate  a characteristic  measure  of  this  time, 
but  no  attempt  was  made  to  model  the  detailed  fluid  dynamics  during  start- 
up . 

Optimum  performance,  in  the  sense  of  minimized  pumping  time,  may 
require  interactive  control  of  the  overall  pumping  rate.  Ideally,  the  flow 
in  the  sender  tank  is  divided  into  two  regions.  A central  gas  core-defined 


3 


by  a single  gas-liquid  interface  is  surrounded  by  a liquid  annulus 
containing  a dispersed  gas  phase.  Sufficient  vortex  rotation  in  the  sender 
tank  must  be  maintained  to  provide  the  radial  pressure  gradient  necessary  to 
transport  the  gas  bubbles  from  the  wall  region  toward  the  interface  and  the 
liquid  from  the  tank  core  region  toward  the  wall.  If  the  rotation  is  too 
weak  gas  bubbles  will  be  drawn  into  the  MIST  unit  along  with  the  liquid.  On 
the  other  hand,  intense  rotation  produced  by  a high  pumping  rate  could 
produce  a sufficiently  high  speed  injection  rate  that  the  injected  fluid 
will  travel  along  the  wall  directly  to  the  MIST  unit  collection  point.  The 
optimum  pumping  rate  can  be  determined  by  modeling  the  flow  in  the  sender 
tank  during  the  continuous  stage  of  operation. 

2 . 1 Early  Transfer  Concept 

During  the  early  stages  of  this  study,  a slightly  different  operational 
mode  was  proposed.  In  this  mode,  the  vortex  motion  was  established  and 
maintained  by  tangential  gas  injection  at  the  wall  as  illustrated 
schematically  in  Figure  1.  Interfacial  momentum  transfer  from  the  gas  phase 
to  the  liquid  phase  was  the  principal  mechanism  sustaining  the  overall 
annular  vortex  flow  pattern  illustrated  in  Figure  2. 

The  transfer  process  was  controlled  by  the  MIST  unit  illustrated  in 
Figure  1.  Two  stages  of  operation  could  be  identified.  During  start-up, 
the  MIST  unit  drew  an  initially  arbitrary  mixture  of  liquid  and  gas  from  the 
sender  tank,  separated  the  phases,  sent  the  liquid  phase  to  the  receiver 
tank,  and  returned  the  gas  phase  to  the  sender  tank  through  the  tangential 
injection  system  (MIST  unit).  The  returning  gases  initiated  rotation  of  the 
tank  contents.  Once  the  annular  liquid  region  shown  in  Figure  2 is 
established,  a valve  closed  in  the  gas  return  line  of  the  MIST  unit  and 
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liquid  was  continuously  transferred  by  the  MIST  unit  to  the  receiver  tank. 
During  this  continuous  stage  of  operation  gas,  displaced  by  the  liquid  in 
the  receiver  tank,  was  transferred  to  the  sender  tank  through  the  tangential 
injection  system  and  sustained  the  vortex  motion. 

After  preliminary  study  established  a general  understanding  of  the 
overall  flow  mechanics,  i.e.,  the  integral  momentum  balances  involved,  the 
conclusion  was  drawn  that  this  early  transfer  concept  should  be  somewhat 
modified.  The  modification  made  was  mainly  in  procedure  rather  than  in 
hardware . 

Conceptually,  the  problem  with  driving  the  liquid  vortex  with  an 
incoming  stream  of  vapor  and  extracting  a stream  of  (mostly)  liquid  is  that 
more  angular  momentum  is  being  withdrawn  from  the  tank  than  is  being 
introduced.  This  is  due  to  the  large  density  difference  between  the  two 
streams.  The  overall  decrease  of  angular  momentum  will  cause  the  vortex  to 
slow  down,  increasing  the  likelihood  that  gas  bubbles  will  be  drawn  into  the 
outlet  stream.  Though  the  flow  in  the  tank  is  complicated,  it  would  seem 
likely  that  the  rotating  flow  could  not  even  be  initiated  using  this 
approach . 

Although  the  incoming  angular  momentum  can  be  increased  by  modifying 
the  inlet  nozzle  geometry,  an  extreme  difference  in  the  nozzle  geometries 
would  be  required  to  maintain  a constant  angular  momentum  due  to  the  large 
density  difference  between  gas  and  liquid  phases.  Additionally,  viscous 
losses  will  further  increase  the  loss  of  angular  momentum.  Loss  of  angular 
momentum  in  the  tank  will  generally  result  in  a decrease  in  angular 
velocity.  Thus,  the  radial  pressure  gradient,  and  therefore  the  centrifugal 
force  on  the  liquid  will  be  reduced. 
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2.2  An  Alternative  Mode  of  Operation 


As  the  result  of  the  general  conclusion  drawn  above,  a somewhat 
modified  mode  of  operation  was  proposed.  Instead  of  the  single  continuous 
process  of  injecting  gas  into  the  tank  while  liquid  is  being  pumped  out,  a 
two-stage  operation  may  be  employed.  The  first  stage  - the  start-up  stage 
(Figure  3) » would  return  all  fluid  - liquid  and  gas  that  is,  pumped  from  the 
outlet  - back  to  the  inlet.  The  greater  density  of  the  injected  fluid 
would  alleviate  the  momentum  loss  discussed  above  and  permit  the  vortex  to 
develop.  Eventually,  vapor  bubbles  would  migrate  to  the  core  of  the  tank 
and  only  liquid  would  be  passed  through  the  pump.  When  the  vortex  has 
sufficient  strength  (or  when  a balance  of  pumping  gains  and  flow  losses 
exists)  and  a certain  phase  separation  has  been  achieved,  the  liquid  can  be 
pumped  to  the  receiver  tank  (Figure  4)  and  vapor  returned  to  the  sender  tank 
(but  to  the  core , not  to  the  perimeter).  In  this  second  stage  - the  pump- 
out stage,  as  the  liquid  is  withdrawn,  the  vortex  will  lose  strength,  but 
since  the  vapor  returns  to  the  core,  there  is  a much  reduced  risk  of 
entraining  bubbles  in  the  outlet  until  the  vortex  strength  is  nearly  spent. 

A conceptual  sketch  of  how  the  two-stage  "start-up"  and  "pump-out" 
model  might  operate  is  shown  in  Figure  5.  Schematically,  the  time 
evolutions  of  the  vortex  strength  and  the  degree  of  phase  separation  of  the 
two-stage  operation  concept  are  plotted  on  an  arbitrary  ordinate  scale  for 
the  purpose  of  describing  the  concept.  The  definitions  and  the  magnitudes 
are  entirely  arbitrary.  Nevertheless,  the  vortex  strength  can  be  thought  of 
as  the  total  angular  moemntum  of  the  fluids  in  the  tank,  whereas  the  degree 
of  phase  separation  can  be  thought  of  as  the  liquid  volume  fraction  near  the 
port  from  which  liquid  will  be  withdrawn.  Other  definitions  can  also  be 
found  in  the  section  on  numerical  results.  Initially,  (i.e.,  at  t = 0)  the 
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system  is  shown  to  have  zero  vortex  strength  and  an  arbitrary  level  of  phase 
separation.  During  the  start-up  stage,  (0  < t < t*)  due  to  fluid  injection 
action,  the  vortex  is  growing  and  its  strength  is  increasing  with  time.  The 
radial  pressure  gradient  produced  by  the  vortex  motion  will  promote  the 
light  gas  to  move  inwardly  and  the  heavy  liquid  to  move  outwardly  and  the 
degree  of  phase  separation  is  thus  increasing  with  time.  As  the  vortex 
motion  reaches  a sufficient  strength,  additional  injection  action  may  no 
longer  be  very  effective  because  higher  velocities  would  be  required.  At 
this  time,  the  injecting  action  can  be  either  reduced  or  stopped.  In  the 
mean  time,  the  separation  process  continues  because  of  the  existing  vortex 
motion.  Once  the  degree  of  separation  has  progressed  beyond  the  critical 
value,  (e.g.,  at  time  = t*) , the  pump-out  stage,  in  which  the  liquid 
transfer  from  the  sender  tank  to  the  receiver  is  done,  can  be  started. 
Critical  separation  levels  (as  one  being  shown  in  Fig.  5)  can  be  determined 
according  to  various  criteria.  As  the  pump-out  (transfer)  stage  proceeds 
(t*  < t < t**)  , the  vortex  strength  will  decrease  by  virtue  of  the  liquid 
being  pumped  out  from  the  annular  layer  at  the  cylinder  wall  and  the  degree 
of  separation  can  either  increase  or  decrease  depending  on  the  strength  of 
the  remaining  vortex  and  the  pumping  conditions.  When  the  degree  of 
separation  reduces  to  the  critical  phase  separation  (e.g.,  time  = t**),  the 
transfer  is  no  longer  practical  due  to  the  small  volume  fraction  of  liquid 
being  able  to  be  transferred,  and  the  transfer  process  stops.  The  liquid 
residual  in  the  sender  tank  when  the  transfer  is  stopped  is  another  design 
parameter  and  should  be  made  as  small  as  possible.  It  can  sometimes  be 
reduced  by  increasing  the  vortex  strength  and/or  by  decreasing  the  pump-out 
rate . 
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As  discussed  above,  a two  stage  operational  mode  is  proposed  in  this 
study.  Figure  6 shows  the  conf iguration  of  the  modified  two-stage  MIST 
Transfer  System.  During  the  start-up  stage,  = 0,  and  . That 

is,  all  fluids  withdrawn  are  returned  to  the  sender  tank.  The  flow  dynamics 
in  the  sender  tank  during  this  stage  are  those  of  the  pertinent  two-phase 
flow  situation.  This  stage  is  employed  to  develop  an  active  vortex  such 
that  the  separation  of  liquid  from  vapor  bubbles  can  be  initiated  and 
ultimately  achieved.  In  the  transfer  stage,  * 0 and  the 

flow  could  be  treated  as  a single  phase-liquid,  pump-out  problem.  To  insure 
an  optimal  liquid  transfer,  a sufficiently  high  degree  of  phase  separation 
during  the  pump-out  stage  needs  to  be  maintained.  For  some  cases,  if  the 
vortex  loses  too  much  strength  either  due  to  a large  pump-out  rate  or  due  to 
friction  during  the  pump-out  stage,  a small  amount  of  the  pumped  liquid, 
could  be  recirculated  back  to  the  sender  tank  to  maintain  sufficient  vortex 
angular  momentum.  Also,  any  gas,  picked  up  and  separated  by  the  MIST 
unit  will  be  sent  back  to  the  core. 

2.3  Start-Up  (Spin-Up)  Stage 

During  the  start-up  period,  since  all  fluid  pumped  out  is  returned  to 
the  tank,  the  mass  of  fluid  inside  the  tank  remains  unchanged.  However,  due 
to  the  action  of  the  external  pump  and  its  injection  action,  the  total 
angular  momentum  (or  vortex  strength)  is  increasing  with  time.  Since  a jet 
at  a larger  radial  location  produces  a larger  angular  momentum  for  the  same 
linear  jet  speed,  the  inlet  nozzle  should  be  located  near  the  perimeter. 
Also,  since  a mixture  with  a higher  liquid  volume  fraction  has  higher  mass 
density  thus  producing  higher  angular  momentum,  the  outlet  nozzle  should 


also  be  located  near  the  perimeter  where  a higher  liquid  volume  fraction 
expected. 

Referring  to  Figure  6,  the  net  change  of  angular  momentum  (or  vortex 

strength)  of  the  system  due  to  injection  action  can  be  estimated  by 

dM 


"tt“  = R(m_V . 
dt  3 j 


m1V) , 


is 


where  R is  the  radius  of  the  tank,  nh  is  the  mass  flow  rate  at  the  port 

location  i , V is  the  vortex  speed,  and  V.  is  the  injected  jet  speed.  Thus, 

as  expected,  in  order  to  increase  the  vortex  strength,  the  jet  speed  should 

be  larger  than  the  vortex  speed,  i.e.,  > V since  m^  = m^.  In  fact, 

unless  the  jet  speed  is  larger  than  the  local  vortex  speed,  the  angular 

momentum  of  the  system  will  normally  decrease  with  time  and  the  desired 

rotational  motion  could  not  even  be  initiated.  This  is  due  to  two  facts. 

First,  during  this  period,  it  is  expected  that  the  lighter  vapor  bubbles 

would  migrate  to  the  core  of  the  tank  and  the  average  density  near  the  wall 

will  increase  with  time.  Thus,  a slight  time  lag  for  the  mixture  to  get  to 

the  inlet  nozzle  from  the  outlet  nozzle  normally  implies  m^  < m^.  Again, 

this  means  more  angular  momentum  is  being  withdrawn  from  the  tank  than  is 

being  introduced  if  V.  = V.  Secondly,  fluid  viscosity  will  further  reduce 

1 

the  total  angular  momentum  of  the  system. 

Some  characteristic  features  to  be  noted  during  the  start-up  stage  can 
be  summarized  below: 

i)  The  fluid  is  a two-phase  one  - where  each  plays  a significant  role; 

ii)  All  of  the  fluid  mixture  which  is  pumped  out  is  returned  to  the 

tank  (Q1  = , Q2  = = 0) ; 

iii)  Both  inlet  and  outlet  nozzles  are  located  near  the  tank  perimeter 
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iv)  The  injected  jet  speed  should  be  larger  than  the  local  vortex  speed 
to  increase  the  net  angular  momentum. 

2. 4 Pump-Out  (Transfer)  Stage 

After  the  vortex  has  sufficient  strength  and  the  spin-up  time  is 
sufficiently  long,  the  separation  of  liquid  and  vapor  can  be  established 
according  to  specified  criteria.  At  this  condition  (i.e.,  t = t*f  see  Fig. 
5),  the  liquid  forms  an  annulus  against  the  wall  of  the  tank  and  around  a 
core  of  vapor,  and  the  second  operational  stage,  the  transfer  stage,  can 
begin.  In  this,  single  phase  liquid,  Q can  be  pumped  out  from  the  tank 
perimeter,  and  vapor  from  the  receiver  tank  returned  to  the  center  core  of 
the  sender  tank,  see  Fig.  6.  During  this  process  the  total  angular  momentum 
will  decrease.  To  maintain  a reasonable  vortex  strength,  a small  amount  of 
high  speed  liquid,  Q , could  be  pumped  in  at  the  perimeter.  Thus,  the 
transfer  stage  is  mainly  a single-phase  flow  with  a central  vapor  core. 

Some  characteristics  of  the  transfer  stage  are: 

i)  There  is  an  interfacial  radius  R^  between  the  vapor  core  and  the 
liquid  annulus;  gas  is  present  only  for  r < Rc;  liquid  is  present  only  for 

r > R ; 

c 

ii)  Inlet  vapor  flows  into  the  center  core  and  outlet  liquid  flow,  Q 
occurs  at  the  tank  perimeter; 

iii)  If  needed,  a small  inlet  liquid  flew,  can  occur  at  the  tank 
perimeter  with  an  inlet  jet  speed  V.  > V; 
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d 2 

iv)  Q„  = 0,  and  transfer  rate  Q0  * tt  (tt  R ) 
4 d dt  c 


2 f ^2  2 

or  R - — dt  + R 

C J TT  CO 


where  R is  the  initial  interfacial  radius  of  the  transfer  stage 
co 

v)  Total  angular  momentum  in  tank  M is 

3 


r Rc  2 r R 

- 2ir[p^  | r V dr  + 

; o R 

c 


r2  V dr] 


and  also 


dM 

— - = o R (Q  V - Q V) 
dt  p2  v3  j 1 


If,  = 0,  and  V = roo,  we  have 


w TTOO  r f \ 

Ma  - - C(Pl  - P2)  R0  - p2R  J 


and 


Thus , we  have 


dMa  ,2  n 

5T  = 'P2R  “Q1 


; _ -2  Qx  [(PX  - P2)  Rc2  - p/] 
“ [(px  * P2)  r]  + PjR4]  if 


An(— ) 
co 

o 


2 V(pl  - p25  Rc  + P2r2] 


J 


(pl  “ P2)Rc  + P2R4 


dt 


1 1 


Since  <<  p^,  the  vortex  velocity  becomes 


J.n(  — ) = — ( 

U)  TT  J 

O 


Q1dt 


R 


dR 


2 2 
0 (R  + R ) 
c 


= -2 


2 2 
R R + R 
CO  c 


0)  R2  + R 2 2 

— = (-p ^f-)  = [1 

“ R2  + R 2 


+ -Jo*L 


Q,dt 


2 2 
ir  (R  + R ) 

CO 


'] 


This  equation  describes  the  decrease  of  the  angular  velocity  of  the  vortex 
as  the  transfer  is  taking  place. 

The  pressure  at  the  tank  wall  is  another  important  quantity.  It  can 

be  estimated  in  the  vortex  flow  where  the  radial  velocity  is  much  smaller 

than  tangential  velocity.  Integration  of  pressure  gradient  equation  ^ = 

2 . 

pvuj  gives 


P 

w 


P 

c 


[plRc' 


)] 


o P2W  ,d2  d 2 . 

P + -r — R - R 

c 2 c 


if  px  <<  pp 


wo  p2  2 
P + — — (R  - R 

c 2 co 


2 f 


Q1dt 


f t 


[1  ♦ 

it  ( R2  + R 


Q. dt  -4 


2) 


co 


where  P^  is  the  pressure  at  center  core.  If  P^  is  kept  constant,  the  wall 
pressure  P^  will  decrease  as  the  transfer  stage  continues.  It  is  also  noted 
that  the  viscous  stress  on  the  cylinder  wall  will  further  decrease  the 
vortex  strength  oj,  and  thus  reduce  further  the  wall  pressure.  Reduced  wall 
pressure  implies  the  pump  will  have  to  pump  harder  to  withdraw  the  liquid. 
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3 . Modeling  of  Start-up  Stage 


The  physical/mathematical  model  for  the  start-up  stage  will  be  based  on 
the  two-phase,  two  fluid  continuum  hypothesis.  The  work  considers  the  flow 
in  the  sender  tank  to  be  specified  by  the  tank  geometry,  fluid  properties, 
nozzle  geometry  and  the  overall  pumping  rate  (the  volume  flow  rate  of  the 
mixture  in  the  tangential  injection  system).  The  flow  in  the  sender  tank  is 
assumed  to  be  two-dimensional  with  no  flow  and  no  gradients  along  the  length 
(axis  of  symmetry)  of  the  tank.  The  fluid  injection  system  is  modeled  as  a 
direct  line  source  for  tangential  injection  and  fluid  withdrawal  is  modeled 
as  a line  sink.  Although  the  cryogenic  fluids  of  most  interest  in  space 
operations  exhibit  significant  interfacial  mass  transfer,  the  fluid  motion 
is  driven  principally  by  the  interfacial  momentum  transfer  from  the  injected 
fluid  mixture  to  the  fluids  in  the  tank.  Therefore,  our  analysis  assumes 
zero  interfacial  mass  transfer  as  would  be  closely  approximated  by  an  air- 
water  system.  Both  phases  are  assumed  incompressible.  Surface  tension 
between  the  fluids  is  assumed  to  be  negligibly  small  compared  to  the  other 
interfacial  forces. 

3.1  General  Formulation 

The  two-fluid  or  separated  flow  model  is  one  of  several  classes  of 
possible  approaches  to  the  analysis  of  two-phase  flow  [2].  It  is  based  on 
separate  conservation  equations  of  mass,  momentum  and  sometimes  energy  in 
two-phase  flow  for  each  phase.  This  model  seems  to  be  especially  suitable 
for  the  uniformly  dispersed  two-phase  flow.  Various  forms  of  the  system  of 
equations  for  two-fluid  two-phase  flow  may  be  found  elsewhere  [ 3 “9 ] . Many 
of  these  are  obtained  via  postulating  macroscopic  qualities  without 
reference  to  microscopic  relations;  others  are  formulated  more  rigorously 
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and  systematically  by  averaging  the  microscopic  equations.  Some  of  these 


are  good  only  for  special  cases  and  are  more  practical;  others  are  more 
general  but  are  usually  too  complicated  to  be  practical.  Nevertheless,  all 
of  them  seem  to  have  retained  the  expressions  for  the  phase  interactions , 
i.e.,  additional  source  terms  or  closure  laws. 

The  basic  time-dependent  fluid  dynamical  equations  consist  of  averaged 
equations  of  motion  for  each  phase,  and  two  continuity  equations  written  in 
terms  of  the  volume  fraction.  Here  we  will  give  briefly  the  formulation  of 
the  system  equations.  For  a detailed  review  of  this  formulation,  the  reader 
is  referred  to  the  paper  of  Drew  [6]. 

In  the  interior  of  each  phase,  the  equations  of  motion  are 


3pk 

* V • P,  V - 0 


(1  ) 


and 


(2) 


At  the  interface,  jump  conditions  are  written  via  conservation  of  mass  and 
momentum  as 


I pk  (Vk  - V.)  . n = 0 


(3) 


k 


and 


Z[^kvk  (vk  - V - V • n - 


(4) 


k 


Here  p is  density 


V is  the  velocity 
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T = - p I + t is  the  total  stress  tensor 


p is  the  pressure 

t is  the  extra  stress  tensor 

I is  the  unit  tensor 

B is  the  body  force  density 

V.  is  the  interfacial  velocity 
a is  the  surface  tension 
k is  the  curvature  of  the  interface 

n is  the  unit  vector  normal  to  the  interface  and  the  subscript  k = 1 
and  2,  denotes  gas  and  liquid  phases,  respectively. 

Although  the  equations  in  the  interior  of  each  phase  are  well  known, 
the  interface  conditions  between  two  phases  are  quite  complicated.  That  is, 
the  essential  difficulty  in  the  two-fluid  formulation  comes  from  the 
presence  of  moving  interfaces  and  the  interfacial  mass  and  momentum 
transfers.  When  one  considers  the  local  (microscopic)  flow  properties,  the 
shape  and  velocity  of  the  interfaces  is  not  known  a priori.  Due  to  the 
complication  of  the  detailed  interfacial  properties,  solving  the  local 
instant  variables  of  the  flow  field  based  on  the  system  of  equations  in  the 
interior  of  each  phase  would  result  in  an  intractable  (and  thus  - hopeless) 
multi-boundary  problem  except  in  the  most  academic  examples.  Thus,  as  in 
the  field  of  turbulence,  the  method  of  the  averaged  process  is  adopted  to 
study  this  complicated  two-phase  flow. 

Referring  to  Ref.  6,  the  averaged  equations  can  be  written  as 

a + = 1 (5) 

1 ^ 
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+ V . a,  V, 
k k 


0 


(6) 


and 


with 


and 


3ak 

at 


dV, 


a.  p.  ,,  = - a.  Vp,  + V . a,  i.  + a,  p,  B.  + M. 

k k dt  k k k k k k k k 


(7) 


M + M2  = 0 


(8) 


dt 


is  the  total  derivative  of  the  velocity  V . 


Here,  all  the  quantities  are  referred  to  as  the  averaged  macroscopic 


values,  a,  is  the  volume  fraction,  t,  is  the  total  effective  extra  stress 
k k 

tensor,  and  is  the  effective  interfacial  force  density  of  the  k-phase.Eq. 
(8)  merely  states  that  the  interfacial  forces  with  opposite  signs  in  the 
momentum  equations  balance  between  the  separate  phases.  That  is,  the 
interfacial  force  is  just  an  internal  force  for  the  mixture  and  will  not 
appear  in  the  equation  of  motion  for  a mixture.  Thus,  from  equation  (7), 
the  mixture  momentum  equation  is 


dVk  = 

z Vk  — = -7p  + 7 • 7 ak  V + z vA  (9) 

k k 

where  the  interfacial  force  has  been  eliminated. 

In  order  to  obtain  the  averaged  system  of  equations  several  assumptions 

have  been  made.  Among  those  are: 

i)  Constant  physical  properties  for  each  phase,  i.e.,  p,  , y,  are 

K K 

constant . 
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ii)  Each  phase  is  treated  as  an  incompressible  fluid  and  there  is  no 
pressure  jump  inside  each  phase. 

iii)  No  phase  change  between  phases,  i.e.,  there  is  no  interfacial  mass 
transfer . 

iv)  Neither  heat  transfer  nor  thermodynamics  is  considered. 

v)  The  physical/mathematical  model  is  based  on  the  two-phase  continuum 
hypothesis.  The  initial  approach  is  to  assume  a bubbly  flow  regime. 

Furthermore,  the  surface  tension  (or  the  contact  pressure)  is 
neglected,  and  there  is  no  pressure  difference  between  the  phases.  Thus, 
only  one  pressure  variable  is  all  we  need,  and  equations  (6)  and  (7),  with 

Px  = P2  - P» 

describe  the  system  dynamics  provided  the  variables  and  can  be 
expressed  in  terms  of  the  local  state  variables. 

Our  averaging  process  eliminates  complex  formulations  for  the 
microscale  fluctuations  around  the  bubbles  or  droplets.  This  process  has 
lumped  the  effects  of  the  interfacial  interactions  between  phases  into  the 

effective  stress  tensor  i,  and  the  effective  interfacial  force,  M,  . In 

k k 

order  to  properly  model  a two-phase  flow,  the  constitutive  equations  for  the 

averaged  stress  tensor  and  the  interfacial  force  density  need  to  be 
given . 

In  general,  both  t and  M could  be  functions  of  various  parameters, 

i\  K 

d?k 

such  as  p , u , V , d^,  , ..  etc.,  in  a very  complicated  way.  In 

modeling  the  averaged  stress  tensor  t , an  analogy  in  certain  aspects 

i\ 
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between  single-phase  and  two-phase  flow  is  assumed.  As  given  in  Refs.  (6, 


10),  the  stress  tensor  due  to  fluctuations,  i is  modeled  by  the  effective 
viscosity  and  the  mean  deformation  tensor  as 

\ ■ 2 (^k + O V 4 (v\  + (10) 


t e t 

where  y , y^  and  y^  = y^  + y^  are  the  dynamic  viscosity,  the  turbulent 
eddy  viscosity,  and  the  total  effective  viscosity,  respectively.  Although 
there  are  several  models  for  the  eddy  viscosity  for  single  phase  flow,  few 
correlations  for  the  eddy  viscosity  in  two-phase  flow  can  be  found  in  the 
literature.  In  reference  10,  the  fluctuations  of  the  turbulent  velocity, 
and  thus  the  turbulent  eddy  viscosity  y t,  in  two-phase  flow  had  been 

r\ 

assumed  to  be  due  to  two  distinct  mechanisms:  (1)  momentum  exchange  of  the 

liquid  phase,  and  (2)  movement  of  the  dispersed  phase.  A proper  modeling  of 

e 

the  total  effective  viscosity  y in  a complex  two-phase  flow  requires 
further  study.  In  this  report  a constant  effective  viscosity  is  assumed. 


The  modeling  of  the  interfacial  force  density  M,  is  much  more  involved 

K 

because  the  interfacial  forces  which  result  from  the  averaging  process  are, 
in  general,  very  complicated.  Many  investigators  [4,  9,  11-15]  have  devoted 
considerable  effort  to  define  or  model  these  interact  ions . The  most  common 
effects  included  are  the  classical  drag  forces,  the  forces  due  to  the 
relative  acceleration  or  the  so  called  added  mass  forces,  the  lift  forces, 
the  Basset  memory  forces,  etc.  Here  we  assume  that  the  interface  forces  are 
dominated  by  the  classical  drag  and  added  mass  terms  as: 


Mi  - 


A.(V,  - V.  ) 
d 2 1 


— (V  - 

dt  ^ 2 


V 


(11) 
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the  first  term  represents  the  classical  drag  forces  and  the  second  term  is 
the  force  due  to  the  relative  acceleration. 

With  the  aid  of  the  constitutive  equations  (or  closure  laws)  for  the 
total  stress  (eq.  10)  and  for  the  interfacial  forces  (eq.  11),  equations 
(5),  (6)  and  (7)  now  describe  the  system  dynamics  of  our  two-phase,  two 
fluid  problem.  The  nine  equations  provide  a set  of  equations  in  the  nine 

variables  (6  variables),  (two  variables),  and  p.  They  are  thus  the 
equations  required  to  provide  a mathematical  description  of  the  two  phase 
flow  system  in  terms  of  local  mean  variables. 

From  equations  (5)  and  (6),  an  important  incompressible  condition  is 
obtained 

Z V . ak  Vk  = V . Z (akV  ) = o (12) 

k k 

Equation  (12)  is  analogous  to  the  important  condition  of  V . V = o for 
incompressible  single  phase  flow.  It  should  be  emphasized  here  that  in  a 

two-phase  flow,  the  condition  of  7 . = o is  no  longer  true  even  though 

the  fluid  is  assumed  to  be  incompressible  for  each  phase  (i.e.,  p = 

K 

constant ) . 

The  equations  of  the  two-phase  model  are  coupled  by  the  interface 
forces.  Although  both  the  drag  and  added  mass  forces  are  providing  the 
coupling  effect  between  the  phases,  the  inclusion  of  the  added  mass  in  the 
interfacial  forces  make  the  momentum  equation  behavior  much  more 
complicated.  This  is  because  the  effects  due  to  these  two  forces  are 
somewhat  different.  In  regard  to  the  system  of  differential  equations,  the 
drag  force  is  just  an  algebra  term  and  does  not  alter  the  characteristic  of 


the  differential  equations;  whereas  the  added  mass  force  involving  the  total 
derivatives  does  modify  the  character istics  of  the  differential  equations. 

By  substituting  the  interfacial  forces,  eq.  (11),  into  equation  (7),  we 

have 


dV!  r 

plal  dt  al7p  + 7 ' alTl  + alpl3l  ” Ad  (V1 


dV  dV, 

V,)  - A (-r-i  - 
2 a dt  dt 


(13) 


and 


dV, 


dV 


p0a 


2 2 dt 


= -a2Vp  + V . a2x2  + a2P2B2  - Ad  ^ - V2) 


A .( 


1 


a dt 


dV; 

FF 


(14) 


where 


d ( ) 
dt 


3(  ) 
at 


+ V . V ( 


is  the  total  derivative  operation. 


Since  the  addition  of  the  added  mass  terms  has  modified  the 
characteristics  of  the  equations,  the  method  of  solving  the  system  can  be 
modified  correspondingly . Before  solving  the  new  system  of  the  differential 
equations,  it  is  logical  to  search  and  rearrange  the  dynamic  equations  (13) 
and  (14)  into  a more  convenient  form.  An  additional  motivation  to  search  a 
new,  more  generalized  form  stems  from  the  important  continuity  condition  of 
the  two-phase  flow,  eq.  (12). 

We  can  rewrite  eqs.  (12)  and  (13)  as: 
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Thus  from  Eq.  (15)  x (p0a0  + A ) + Eq.  (16)  x A , we  have 

c.  d.  a a 

d?i 

Cplp2ala2  + Aa  CplV  p2°2)]  dt~  ’ “ (al°2p2  + Aa)Vp 
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and  from  Eq.  (15)  x A + Eq.  (16)  x (p_a,  + A ) we  have 
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Equations  (17)  and  (18)  can  be  written  into  a generalized  equation  as 
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where  the  generalized  coefficients  are 


Cpl  = (al°2p2  + Aa)/<P  > 
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Thus,  equations  (5),  (6)  and  (19)  are  the  generalized  system  of  equations. 

The  major  advantages  of  writing  the  momentum  equations  in  the 
generalized  form  of  eq.  (19)  are  the  following.  First,  the  generalized 
equation  is  similar  to  the  well  known  Navier-Stokes  equation.  Second,  by 
using  eq.  (19)  and  the  continuity  condition  for  two-phase  flow,  eq.  (12),  a 
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Poisson  type  equation  for  pressure  field  can  be  obtained.  Thus,  the 
existing  techniques  can  be  used  to  solve  this  system  of  equations,  and  the 
pressure  field  can  be  determined  through  the  local  independent  variables. 
To  do  this,  taking  the  divergence  of  the  sum  of  equations  (19),  we  have 
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(ctlCpl 
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since 

3 a.  V 

V . Z — - = o from  eq.  (12). 


3.2  Interfacial  Force  Coefficients 

The  modeling  of  the  interfacial  forces  is  one  of  the  most  challenging 
tasks  facing  those  involved  with  the  two-phase  flow  problems.  As  discussed 
before,  the  classical  drag  and  added  mass  terms  are  the  most  common  and 
important  contributions  cited  in  the  literature.  Our  discussion  which 
follows  will  be  limited  only  to  these  two  forces.  In  fact,  even  for  these 
classical  forces,  there  are  various  expressions  used  in  the  literature  [4, 
6,  9,  12,  14].  However,  all  of  these  expressions  are  developed  and  tested 
only  for  dilute,  two-phase  cases. 

Before  modeling  the  interfacial  forces,  it  is  important  to  discuss  the 
average  bubble  or  droplet  size.  The  reason  for  this  is  that  the  two-phase, 
two  fluid  model  basically  assumes  the  flow  field  to  be  in  the  uniformly 
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bubbly  or  in  the  uniformly  dispersed  droplet  flow  regime  and  that  the 
interfacial  forces  are  greatly  dependent  on  the  size  of  bubbles  or  droplets. 
The  importance  of  the  bubble  sizes  on  the  drag  forces  can  easily  be  seen  in 
the  following  discussion.  To  be  more  realistic,  the  average  size,  d^  should 
be  allowed  to  vary  in  the  flow  field  and  in  general  should  be  functions  of 
various  parameters  including  the  local  turbulence  intensity,  mixing,  and  the 
volume  fraction,  a , since  all  these  quantities  will  affect  the  break-up  and 

K. 

coalesce  of  the  bubbles  or  droplets.  To  discuss  the  volume  fraction  effect, 
where  n,  is  the  particle  (bubble)  number  density  of  phase  k,  we  have  a,  = n, 

IK  K K 

3 Y 1-3Y 

Trd.  /6.  And  if  d,  - a,  (i.e.  d,  = d , a,  ) then  n,  - a,  . Since  one 
k kk  kokk  kk 

might  expect  that  both  n.  and  d,  will  increase  as  a,  increases,  we  have  0 < 

k k k 

Y < 1/3.  The  size  variations  as  functions  of  gas  volume  fraction,  a,  for  Y 
= 0.2  and  d^  = dQ2  = 1 cm  are  shown  in  Fig.  7.  Although  these  are  the 
typical  values  used  in  the  test  runs,  it  is  noted  that  other  expressions  can 
easily  be  put  into  the  program.  After  having  the  size  information,  the 
interfacial  forces  can  be  modeled.  For  dilute  bubbles  dispersed  in  a 
continuous  liquid,  the  drag  and  added  mass  coefficients,  A^  and  A^  are 
generally  expressed  as  [6,  12,  14] 
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Here  the  bubble  particles  are  assumed  to  remain  spherical, 

CQ  = 24  (1  + 0.1  Re374  + . . . )/Re 


24 


via  the  Oseen  drag  coefficient  [11],  where 


Re  = 


padi I vi~v2 I 


2m 


, and  u2m  - u2/a2. 


Here,  d^  is  the  averaged  diameter  of  the  dispersed  bubbles.  If  the  Stoke ’s 
drag  is  used  (i.e.  = 24/Re)  then  the  drag  coefficient  becomes 


1 ^y2al 

'di  = 72 

dl  a2 
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which  is  independent  of  the  velocity  field.  Also,  since  the  upper  bound  of 
the  added  mass  coefficient  should  not  exceed  a2p2  which  corresponds  to 
carrying  all  of  the  continuous  liquid  with  the  bubble  during  its 
acceleration,  a modified  added  mass  coefficient  is  assumed  to  be 


/I  3 \ cup- 

alp2  2 + 2 al  * 2 


(24) 


On  the  other  hand,  if  one  considers  that  the  two-phase  flow  consists  of 
liquid  drops  dispersed  in  a continuous  gas  flow  field,  a second  set  of  the 
coefficients  can  be  expressed: 
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where  the  maximum  percentage  of  0.8  for  droplets  in  gas  is  assumed. 

Now  we  have  two  sets  of  expressions  for  the  added  mass  and  drag 

coefficients,  but  no  existing  criteria  to  determine  which  set  of  equations 

should  be  used  for  a general  two-phase  flow  problem.  As  a matter  of  fact, 

neither  set  of  the  expressions  is  expected  to  be  good  for  the  entire  range 

of  volume  fraction.  For  a dilute  gas  dispersed  in  a liquid  (i.e.  0), 

the  coefficients  should  approach  those  of  the  first  set  (i.e.  A,n  and  A ,) 

dl  al 

while  for  a dilute  liquid  dispersed  in  gas  (i.e.  a ->  1)  the  coefficients 

should  approach  those  of  the  second  set  (Ajn  and  A _).  In  general,  the 

d2  a2 

value  of  could  be  any  value  between  0 and  1.  Here  we  will  model  the 
coefficients  by  a weighting  technique,  and  express  the  coefficient  as 
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Adl  Wdl  + Ad2  Wd2 
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where  w , , w w._  and  w.*  are  the  weighting  functions  to  be  selected, 
al  a2  dl  d2 

Normally,  we  set  w _ + w _ = w.,  + w_,^  = 1 and  all  the  weighting  functions 

J al  a2  dl  d2 

are  functions  of  the  volume  fraction  ct^.  In  the  model  here,  we  choose 
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with  n = n , = 4. 
a d 

Fig.  8 shows  the  dilute  drag  coefficients  A^  and  A^  and  the  modeling 
drag  coefficients  A^.  Fig.  9 shows  the  corresponding  added  mass 
coefficients.  The  modeling  coefficients  approach  those  for  dilute 
conditions  as  should  be  the  case.  It  should  be  noted  again  that  different 
models  for  the  coefficients  can  be  easily  adopted. 


3.3  Equations  in  Polar  Coordinates 

Due  to  the  nature  of  the  problem,  it  is  more  convenient  to  write  the 
system  of  equations  in  polar  coordinates.  For  a two-dimensional  flow  in 
polar  coordinates,  some  of  the  useful  relations  are: 


V(  ) _ 21i  F + 2L!  9 
3r  r36 


- 1 3r(  >r  3(  )9 

V . < ) - - - ' * 

r 9r  r98 


V . 


r 1 ^ Vr  3(  V ( )ea  - 
(-  , [I  . — + —Jr:  — _®i]  r 


r 9r 


r99 


' ) 


+ [l  3r)  V9  ♦ 3^99  - ? 


'r  9r 


r30 


dV,  9a,  V, 

k k k _ 77  77 

ak  dt  " 9t  V ‘ akVkVk 


(31) 


27 


3V 


rrk 


= 2 y, 


rk 


3r 


Tr0k 


/r) 


1 0 6k 


2 \ Vrk/r 


(32) 


Thus,  equations  (6)  and  (19)  become 
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for  k = 1 and  2. 
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For  the  axisymmetric  case,  — = 0,  and  we  have 
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for  k = 1 and  2.  The  incompressible  condition  of  eq.  (12)  reduces  to  a^V  ^ 

+ a2Vr2  = Qr/r»  where  Qr  is  the  net  radial  volumetric  flow.  In  the  start-up 
stage,  = 0. 


3.  Initial  and  Boundary  Conditions 

Eqs.  (33)*(35)  are  the  governing  differential  equations  for  the  two- 
phase  fluid  motion  in  cylindrical  coordinates  under  zero  gravity  conditions. 
To  completely  define  the  system,  we  need  to  specify  both  the  initial  values 
and  the  boundary  conditions. 
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The  program  was  developed  to  accept  arbitrary  initial  conditions.  In 

addition,  a few  classical  initial  velocity  fields  are  also  provided.  All 

these  classical  velocity  fields  have  zero  radial  velocities,  i.e.  V =0  and 

r 

their  tangential  velocities  are: 

1)  Pure  rotation 
V0  = ru, 

where  oj  is  the  circular  rotation  speed. 

2)  Hammel-Oseen  Vortex  (H.O.) 


ve  - 


1 .398  V , -2, 

p ^ _ g-l. 25643  r ) 


3)  G.  I.  Taylor  Vortex  (G.I.T.) 


,,  ,,  - (l-r2)/2 

V = V re 
0 p 


In  H.O.  and  G.I.T.  vorticies,  r = r/r  , r is  the  vortex  core  or  the 

P P 

location  of  the  maximum  velocity,  and  V^  is  the  maximum  velocity.  Also  when 
the  classical  vortex  is  used,  the  velocities  for  both  phases  are  assumed  to 
be  the  same. 

The  specification  of  the  initial  condition  is  rather  straightf orward , 
while  the  boundary  conditions  are  much  more  involved.  In  general,  the 
mixture  could  be  pumped  in  and  out  of  the  system  in  a very  complicated  way. 
If  determining  the  detailed  flow  field  is  the  primary  interest,  detailed 
boundary  conditions  should  be  included.  Here,  we  have  tried  to  simplify  the 
situation  somewhat  and  model  the  boundary  condition  via  global  effects  since 
only  the  global  dynamical  behavior  of  the  two-phase  flow  in  the  tank  is 
expected  to  be  useful  from  the  model.  For  example,  at  the  boundary  location 
(i.e.  at  the  inlet  or  outlet  nozzles)  the  pumped  flow  conditions  are 


specified  by  the  average  volume  fraction  of  the  fluid  mixture  a , the 
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average  pump  volume  density  q.,  the  average  jet  speed  V.,  and  the  effective 

r J J 

pump  region,  PL  = j hdR  = hAR,  where  AR  is  the  nozzle  width,  and  h is  the 
effective  jet  penetration  depth.  For  axisymmetric  cases,  h = 2mr.  The 


global  values  of  V.  and  q.  could  be  estimated  from  pump  flow  rate  Q.  and  the 

1 

nozzle  geometry  as:  V.  * Q./AR,  and  q.  = Q./H.  = V./h.  Thus,  the  volume 

11  1111 

source  rate  for  the  k-phase  at  the  nozzle  location  is 


% = ak  q.,  («0) 

and  the  momentum  flux  density  for  the  k-phase  at  the  jet  location  is 

“k = = pkVjY  (41) 

where  n is  the  unit  vector  of  the  jet. 

In  general,  these  source  terms  should  be  added  to  the  system  of 
equations  at  the  nozzle  location.  That  is,  a volume  source  term  of  rq  for 

K 

the  continuity  equation  (33)  and  a generalized  momentum  flux  density  of 
rct^  I ft^n  for  eqs.  (3*0  and  (35)  should  be  added  at  the  nozzle  location. 
Since  the  cylindrical  boundary  is  fixed,  the  net  inlet  flow  should  be 

zero , 

Eq.H.  = I Q,  = 0 (42) 

k j j 

The  summation  is  taken  over  all  the  inlet  and  outlet  nozzles. 

In  reality,  a source  or  sink  is  arranged  at  the  nozzle  location. 
However,  during  the  start-up  stage  the  mixture  pumped  out  is  injected 
immediately  back  into  the  tank  at  the  nearby  location.  The  net  volume  or 
mass  in  the  system  is  effectively  unchanged  except  for  the  net  change  on  the 
angular  momentum.  Thus,  the  pump  system  (withdrawal  and  injection)  during 
the  start-up  stage  acts  as  a body  force  on  the  mixture  at  the  nozzle 
location.  The  net  momentum  gain  during  the  start-up  stage  is  thus  the 
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momentum  introduced  into  the  system  minus  the  local  momentum  pumped  out. 
Therefore,  we  will  model  this  pumping  dynamic  by  body  forces  without 
considering  the  mass  transfer.  That  is,  in  equations  (34)  and  (35), the  body 


force  density  a p B will  be  replaced  by  the  net  momentum  gain,  Q -a,  p q.V 

K K K K K K J 

and 


VkV  = w/(Y ' v)  = J^A  CY  - v) 


(43) 


will  be  assumed  at  the  nozzle  location.  So  the  total  net  momentum  gain  due 
to  the  pumping  is 

E E | a.p.V.  (V.n  - V)  dr  = E E a,p.Q.  (V.n  - V)  (44) 
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Besides  the  pump  flow  condition,  the  boundary  condition  at  the  tank  wall  is 
also  required.  In  general,  a non-slip  condition  at  the  cylindrical  wall 


should  be  used,  (i.e.  = 0 at  r = R).  If  the  detailed  flow  field  is  not 

expected  to  greatly  affect  the  global  phenomenon  of  the  fuel  transfer  of  the 
tank  system,  the  shearing  stress  at  the  wall  could  be  simulated  by  a 
momentum  sink  (a  negative  source)  near  the  wall.  For  example,  these  wall 
stresses  could  be  modeled  as 
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where  t , is  similar  to  a momentum  sink  at  the  wall,  since  i , < 0.  The 

wk  wk 

condition  of  no  radial  velocity  at  the  boundary  (i.e.  V = 0 at  r = R)  will 
be  maintained.  To  avoid  a numerical  singularity  at  the  origin,  a finite 
center  core  region  was  excluded  from  the  analysis.  At  this  inner  boundary, 
a free  surface  with  no  radial  velocity  will  be  assumed  at  all  times. 
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3 . 5 Numerical  Method 


Equations  (5),  (33),  (3*0  and  (35)  comprise  the  system  of  equations  for 
an  axisymmetric , two-phase  fluid  flow.  The  complete  solution  of  these 
complicated  equations  can  only  be  obtained  through  numerical  methods.  An 
improved  Lax-Wendrof f , two-step  scheme,  (also  referred  to  as  MacCormack’s 
method)  [16,  17,  18]  is  adopted  for  solving  this  time-dependent  problem. 

This  non-centered  differencing  scheme,  using  a full  step  backward  prediction 
and  forward  correction  version,  requires  no  explicit  artificial  viscosity  if 
a proper  stability  condition  is  satisfied.  Using  this  technique  for  solving 
fluid  flow  problems  is  very  efficient  and  has  been  in  widespread  and 
successful  use  for  some  time.  It  is  good  both  for  the  time-accurate 
computation  of  steady  and  unsteady  flow  problems.  The  general  features  of 
the  scheme  are:  i)  its  explicitly  conservative  form,  ii)  it  is  a two-step 
predictor-correction  type,  iii)  it  is  three  point,  two  level  - that  is,  the 
solution  of  f^+1  at  level  n+1  depends  only  on  three  values  of  f*?  at  level  n, 
and  iv)  it  is  second-order  accurate  in  time  and  in  space. 

3.5.1  Equations  Utilized 

Since  there  are  five  independent  variables  (e.g.,  V , vr2»  Vqx»  Vq2 
and  a1  ) for  the  system,  only  five  independent  equations  are  needed  to  define 
the  time  dependent  problem.  The  equations  selected  are  one  continuity 
equation  and  four  momentum  equations.  The  final  solution  is  thus  to  find 
the  correct  pressure  field  so  that  the  other  continuity  equation  is 
satisfied.  Also,  the  condition  of  equation  (5),  = 1 is  used  all  the 

time . 

For  using  the  MacCormack’s  numerical  technique,  the  system  of  equations 
is  expressed  in  the  conservative  form  as: 
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Here  the  subscripts  (t  and  r) 
t and  r,  respectively,  and  W, 
elements  defined  as  follows. 
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denote  partial  differentiation  with  respect  to 
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The  components  of  F: 


The  pressure  components: 
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The  components  of  gG^ : 
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al  [ V91  + Cdir  (Vrl  ‘ Vr2}  + ^ C1Z  (0LZplBr<lr  “ a2,T96<l)] 


a2  ^ V92  + Cd2r  (Vr  1 ~ Vr2)  + ^ C2£  (ailPZBr<lr  “ a2,T992,^ 


al  [“Vrl  V91  + Cdir  (V91  “ Vd2)  + 1 C1Z  (a2,p£3 Blr  ~ 

&=1 


a2  ^ Vr2  V92  + Cd2r  (V01  “ V92)  + ^ C22,  (alPlB 9V  ~ Vrel^ 


(51) 


All  the  components  of  F,  P , 
the  components  of  W which  are  the 


gG^  and  S can  be  regarded  as  functions  of 
independent  variables. 
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3.5.2  Numerical  Scheme 

One  way  to  solve  equation  (46)  numerically  is  to  solve  the  pressure 
field  first  from  equation  (22)  then  integrate  equation  (46)  to  obtain  the 
solution  at  the  new  time  step.  To  ensure  an  accurate  result,  the  solution 
so  obtained  should  be  used  as  a check  to  satisfy  both  the  continuity 
equations  in  addition  to  the  phase  definition  of  a,  + = 1.  In  general, 

the  continuity  equations  probably  will  not  be  satisfied  by  the  solution  due 
to  the  numerical  approximation.  Thus,  an  iteration  on  the  pressure  field  or 
other  treatment  would  be  needed  to  attain  a desirable  accuracy. 

Another  numerical  procedure  that  has  been  used  very  extensively  [19] 

will  be  adopted  here  for  our  purposes.  Basically,  the  procedure  first 

solves  the  equations  with  no  pressure  terms,  then  adds  the  pressure 

correction  back  to  the  solution.  Thus,  from  equation  (46),  we  could  have 

,n 


- W „n  . _^n  , nn 

— — = F + gG  + S 

At  r ° r 


and 


n+l 

W - W* 


= P 


52) 


(53) 


At  r 

where  W*  is  the  intermediate  solution  of  W without  pressure  terms  at  the  new 
time  step  n+l,  the  superscript  n and  n+l  denote  the  values  at  the  initial 
time  and  the  new  time  steps,  respectively,  and  the  elements  of  have  the 
similar  format  as  those  of  equation  (49),  i.e., 
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t 


at.  c _ 

1 pi 


a2Cp2 


0 

0 


(54) 


The  overhead  symbol  ~ here  indicates  the  corrected  terms  of  the  pressure 
field.  This  is  used  to  distinguish  the  pressure  from  those  of  the  regular 
time  steps  n or  n+1. 

From  the  continuity  condition  of  equation  (12),  we  have 

W + W n_  W + Wa  n__ 

V . [(  - ■-■-■)  r + ( -)  9]  = 0 (55) 

r r 

for  each  time  step.  Substituting  eqs.  (54)  and  (55)  into  equation  (53).  We 
have  the  Poisson  type  equation  for  pressure  field  as 

* * * * 

W + W W + W 

V . [(a.C  . + a,C  ,)  Vp]  it  = V . — -)  8]  (56) 

1 pl  2 p2  r r 

Thus,  the  numerical  procedure  is:  1)  to  get  W*  by  solving  the  system  of 
equations  (52)  without  the  pressure  gradient  terms;  2)  to  get  pressure 
gradient  Vp  by  solving  eq.  (56);  and  3)  then  to  correct  the  solution  by 
using  equation  (53). 
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For  axisymmetric  cases,  equation  (56)  reduces  to 


iP 

3r 


* * 
W + w2 


(alCpl  + a2Cp2}  rAt 


and 


= 0 
30 

And  the  solution  after  the  pressure  correction  is 


w n+1  = -w  n+1 
1 2 


a0C  W,  - a C W 
2 p2  i 1 pi  2 

a C + a0C  _ 

1 pi  2 p2 


(57) 


(58) 


n+1  * n+1  * n+1  * 

w3  = w3,  w4  1 = w4,  w5  = w5. 


The  fundamental  theory  of  the  MacCormack’s  scheme  is 
below. 


For  second  order  accuracy,  the  solution  of  eq.  (46) 


. .1  iro  ..  tT  o (At)  Tr  o 
W = W + At  W^.  + — - — W 


briefly  given 


could  be  written  as 


where 


and 


TTo  At  ..  o At  /Tr  o . . TT  oN 

= W + r Wt  - — (Wt  + At  Wtt  ) 

= | (W°  + At  W °)  + | (W°  + At  W P) 

= | (WP  + W°)  (59) 

WP  = W°  + At  W^0  is  the  predicted  value, 

W°  = W°  + At  WtP  is  the  corrected  value. 
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The  superscripts  denote  the  time-level  of  the  information  and  subscripts 
denote  the  partial  derivative  with  respect  to  either  time  t or  space  r. 
Specifically,  superscripts  0 and  1 are  the  initial  and  the  completely 
advanced  time  (here  two  steps)  plane;  p and  c are  the  predicted  (1st  step) 
and  corrected  (2nd  step)  time  plane.  Thus,  W^0  is  the  time  derivative  of  W 
evaluated  at  the  initial  time,  and  is  time  derivative  of  W evaluated  at 
the  predicted  time. 

Fig.  10  shows  the  diagram  of  the  two  step  difference  scheme  used  in  the 
computer  program.  Due  to  the  difference  scheme,  the  spatial  location  after 
each  step  in  time  is  a half  grid  off  from  the  original  one.  Thus,  the 
spatial  offset  which  resulted  from  a backward  predicting  step  will  cancel 
with  those  of  the  forward  correcting  step. 

Numerically,  the  predicted  values  are 


uP  - — (TjP 
l 2 1-1/2 


WP  ) 

Wi+l/2 


(60) 


where 


Wi-S/2  “ I (Wi-l  + «i0)  + At  Ht° 


f O Ox 

= I (w  0 ♦ W °)  + ^ [(F  ° - F ° ) + (G  ° - G ° )] 

2 ^ l-l  i ; Ar  l i-l;  2 k i i-l;j 


+ f {Sl  + Si-P  + At  Pi-l/2 


(61) 


and  the  corrected  value  is  evaluated  at  the  predicted  time  place,  that  is  at 


wm/2-  Thus 
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Wc  = W.°  + At  W 
l l 


P 


t 


2 


+ 


) + At  P.C 


(62) 


p c 

Here  P.  , and  P.  are  the  pressure  correction  terms  as  discussed 

1-1/2  l 


before  (eq.  52-58).  Thus,  for  each  time  step,  the  advance  is  carried  out  in 
two  steps:  a full  step  backward  predictor,  eq.  ( 6 1 ) and  then  a forward 
corrector,  eq.  (62).  As  indicated  in  the  diagram,  the  subscript  i in  the 
regular  mesh  spatial  location  at  which  solution  is  to  be  advanced,  i + 1 is 
the  spatial  location  of  regular  mesh  points  immediately  to  the  right  and 
left  of  the  location  i,  i +_  1/2  is  the  location  midway  between  i and  i + 1 
or  between  i - 1 and  i at  the  predictor  plane.  Thus,  for  each  time  step  as 
the  procedure  advanced,  the  outermost  data  points  at  the  boundary  are  not 
updated  through  the  numerical  scheme.  The  values  at  the  boundary  are  to  be 
given  through  some  suitable  boundary  conditions.  The  numerical  procedure 
utilizes  a uniformly  preselected  spatial  mesh  and  variable  time  increment. 

To  avoid  a singularity  at  the  center  of  the  core  region,  a finite  radius  R. 
is  used  for  the  inner  boundary.  The  tank  radius  R is  the  outer  boundary. 

The  time  step  is  determined  at  each  time  step  to  ensure  numerical  stability 
[20].  For  a finite  grid  size  Ar , the  maximum  time  step  At  is  given  by 


At 


k 


Ar 
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where  k = 1 and  2.  The  minimum  At^  (with  some  rounding  off)  is  used  for  the 
time  step.  Normally,  the  technique  with  the  time  step  condition  of  equation 
(63)  gives  fairly  good  numerical  stability.  However,  in  critical  conditions 
numerical  damping  can  be  added  either  for  damping  oscillations  due  to  large 
gradients  or  for  accelerating  the  calculation  by  increasing  the  time  step. 

A damping  factor,  D thus  was  added  in  the  program  as 

W . 1D  = W.1  (1-D)  + (W1  + W1  - W1)  D (64) 

l l l-l  l+l  l 

where  is  the  value  obtained  based  on  the  two-step  scheme,  and  W1D  is  the 
value  after  the  damping  factor  D is  added.  A typical  value  of  D is  0.2. 

All  the  examples  presented  in  the  following  section  are  based  on  D=0.  That 
is,  no  damping  factor  is  used. 

The  computer  program  was  written  in  a Fortran  77  based  computer  code. 
The  code  will  permit  evaluation  of  the  effects  of  various  parameters  which 
control  the  fluid  dynamical  behavior.  These  include  tank  size,  fluid 
properties,  such  as  density  and  viscosity,  etc.,  characteristic  gas  bubble 
and  liquid  drop  sizes,  and  relative  location  of  injection  nozzles.  Also, 
the  program  was  developed  and  arranged  in  such  way  that  different  modelings 
can  be  easily  adapted. 

3 . 6 Numerical  Results 

The  flexibility  and  ease  of  changing  test  conditions  are  the  major 
advantages  of  a numerical  test  as  opposed  to  an  experimental  test.  The 
computer  program  is  developed  to  accept  various  conditions  and  model 
parameters.  Fluid  properties,  (p,  , u,  ) , initial  volume  fraction 

K K 

distributions  (a^),  initial  velocity  distributions,  boundary  and  injection 
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conditions,  Reynolds  number  and  several  modeling  constants  are  among  the 
input  quantities.  Various  cases  have  been  tested  during  the  course  of  the 
program  development.  To  limit  the  size  of  this  report,  only  a few  typical 
cases  will  be  presented.  The  common  conditions  for  all  the  test  runs 
presented  here  are  as  follows: 

1)  Uniformly  distributed  initial  volume  fraction  (initially,  gas  volume 
fraction  = 0.25  everywhere). 

-R  -?  R 

2)  Air/water  case:  p /p^  = x 10  , y /y^  = 1-13  x 10  , = 10 

R -R 

kg/m  , and  = 1.514  x 10  kg/m-s. 

0 2 

3)  Bubble/droplet  size  obeys  the  power  law  of  d,  = 0.01  x a, 

4)  Jet  opening:  0.85  < r^  < 0.95 

5)  Weighting  powers  for  drag  and  added  mass  coefficients:  n^  = n^  = 4 

5 

6)  Character istic  Reynolds  number,  Re  = p V R/y£  = 10 

7)  Length  scale  = tank  radius  R = 1 meter,  Velocity  scale  = 0.1514  m/s,  Time 
scale  = 6.605  s,  and  Pressure  scale  = 22.922  P . 

a 

The  other  test  conditions  are  summarized  in  Table  I.  These  test  runs  assume 

either  an  initially  pure  rotational  flow  (linear  swirl  velocity 

distribution)  or  that  the  fluid  is  initially  at  rest.  The  initial  angular 

rotation  speed  is  w (radians  per  unit  time).  The  averaged  mean  injection 

speed  across  the  jet  opening  is  V^ . To  specify  the  boundary  condition  at 

the  tank  wall,  I is  used  - with  I =1  being  a non-slip  condition  and  I 
w w w 

0 being  a free  surface  condition  with  no  shearing  forces.  In  either  case, 
there  is  no  radial  velocity  at  the  wall  boundary.  In  addition,  I is  also 
used  to  modify  the  initial  velocity  condition  near  the  wall.  For  1=1,  a 
factor  of  (1  - r)^*^  was  included  on  the  initial  tangential  velocity  field 
when  a simple  vortex  flow  is  assumed.  This  is  assumed  to  simulate  the  power 
law  velocity  distribution  near  the  wall.  To  avoid  the  effects  of  numerical 
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Run 


Speed,  V. 

J 


I 


w 


Vui 


VU2 


No . 


148 

0)  = 1 

1 

0 

115 

1000 

149 

i — t 

it 

3 

1 

0 

1000 

1000 

150 

0)  = 1 

0 

0 

115 

1000 

151 

at  rest 

1 

10,  (t  < 3) 

115 

1000 

0,  (t  > 3) 

152 

at  rest 

1 

10,  (t  > 0) 

115 

1000 

153 

at  rest 

1 

10,  (t  > 0) 

1000 

1000 

154 

at  rest 

1 

5,  (t  > 0) 

115 

1000 

For  all  the  above  cases  the  following  conditions  prevail:  1)  Uniformly- 
distributed  initial  gas  volume  fraction  of  = 0.25;  2)  p /p^  = 1*293  x 

-2-2  02 
10  , y /y  = 1*13  x 10  ; 3)  Bubble/droplet  size  obeys  d^  = 0.001  x * ; 

4)  Jet  opening:  0.85  < r.  < 0.95;  5)  n = n = 4;  6)  Re  = 10^;  7)  Air-water 

J a q 

3 3 -3 

fluids  (p2  = 10  kg/m  , y^  = 1.514  x 10  kg/m-s),  Length  scale  = tank 
radius  = 1 meter,  velocity  scale  Vg  = 0.1514  m/s,  Time  scale  = 6.605  s,  and 
pressure  scale  = 22.922  P . 

Table  I.  Conditions  of  Test  Runs. 


singularities  at  the  origin,  a finite  center  core  region  was  excluded  from 
the  analysis.  For  the  data  presented  here  the  inner  radius  is  10%  of  the 
tank  radius.  At  this  inner  boundary,  the  free  surface  with  no  shearing 
forces  is  assumed  at  all  times.  The  eddy  viscosity  factor  is  y^/y^.  For 

liquid  phase,  y^/y^  = 1000  was  assumed  for  all  cases,  while  two  values  of 

t t 

y /y  were  tested  for  gas  phases.  The  value  of  y /y^  = 115  corresponds  to 

the  condition  of  the  same  kinematic  eddy  viscosity  for  both  phases  (i.e.,  v? 

= v^).  This  is  based  on  the  expression  of  the  eddy  stress  that  pvv  = ytVV. 

As  discussed  by  Drew  [6],  if  both  phases  have  similar  velocity  fields,  the 

kinematic  eddy  viscosities  (v^  = y,  / p , ) representing  both  velocity  scales 

and  mixing  lengths  of  large  eddy  processes  should  be  approximately  equal. 

If  this  is  the  case,  then  the  eddy  viscosity  factor  y^/y  = . 

For  air  water  case,  if  y^/y^  = 1000,  then  y^/y^  * H5. 

Thus,  as  shown  on  Table  I,  test  runs  148  to  150  are  for  an  initially 
pure  rotational  flow  with  no  flow  injection.  Test  runs  151  to  154  begin 
with  the  fluid  at  rest  initially  but  with  flow  injection.  Although  the 
initially  pure  rotational  flow  cannot  exist  in  reality,  it  is  a good  way  to 
test  the  essential  physical  mechanisms  of  vortex-induced,  phase-separation. 
On  the  other  hand,  the  vortex  flow  induced  by  the  flow  injection  of  tests 
151  to  154  are  good  simulations  to  the  actual  situation.  These  spin-up 
tests  are  thus  good  examples  to  study  the  fluid  dynamical  behavior  of  the 
two-phase  fluids  in  the  sender  tank. 

Figs.  11  to  15  show  the  velocity  vector  distributions  for  various  test 
runs.  The  test  results  for  #149  and  #154  are  not  shown  here  because  the 
result  of  #149  is  practically  identical  to  that  of  #148,  and  that  of  #154  is 
very  similar  to  that  of  #152.  In  these  figures,  Fig.  (a)  shows  the  gas 
phase  and  Fig.  (b)  shows  the  liquid  phase.  At  this  point,  it  should  be 
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noted  that  the  actual  flow  fields  are  all  axisymmetric  and  the  angular 
coordinate  replaced  by  the  temporal  variation  in  the  displays  is  done  only 
for  ease  of  observation.  To  produce  these  displays,  the  radial 
distributions  at  the  different  times  are  oriented  along  equally-spaced 
angles.  This  angular  separation  arrangement  should  not  be  interpreted  as  an 
angular  variation  in  the  flow  field  for  the  results. 

In  Figs.  (13),  (14)  and  (15),  the  region  of  injection  is  shown  in  the 
annular  region  between  two  dashed  circles. 

These  velocity  vector  fields  indicate  all  flows  are  primary  in  angular 
rotation  with  gas  phase  tending  to  move  inward  and  liquid  phase  trying  to 
move  outward,  as  expected.  Since  the  wall  shear  stress  will  slow  the  flow 
movement,  the  vortex  flow  will  diminish  eventually  if  there  is  no  injection 
introduced,  as  Figs.  (11)  and  (13)  indicate.  The  effect  of  eddy  viscosity 
factor  can  also  be  seen  from  Figs.  (14)  and  (15).  Fig.  15  shows  the 
velocity  distribution  is  more  uniform  and  extends  more  into  the  center 
region  than  that  of  Fig.  14.  This  is  due  to  the  larger  gas  eddy  viscosity 
factor,  y^/y  on  #153  than  on  #152.  Since  the  main  mechanism  to  establish  a 
flow  rotation  outside  the  region  of  injection  is  the  shearing  stress 
distribution  that  results  from  the  finite  effective  viscosity  flow  having 
larger  eddy  viscosity  should  transfer  and  diffuse  the  locally  injected 
momentum  more  effectively.  Although  the  effect  of  eddy  viscosity  can  be 
found  in  the  injection  cases,  it  shows  practically  no  difference  between 
cases  of  #148  and  #149.  This  is  because  in  these  cases,  the  flows  are 
already  in  pure  rotation  and  thus  the  transfer  or  diffusion  of  the  momentum 
is  not  very  important. 

Figs.  16(a) , (b ), (c)  and  (d)  show  more  details  of  the  velocity 
distributions  of  the  test  #148.  To  keep  the  figure  clear,  most  of  the 
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following  figures  will  show  only  the  data  at  six  times,  t = 0,1, 2, 3,4,  and 
5.  The  tangential  velocities  shown  in  Figs.  16(a)  and  (b)  for  gas  and 
liquid  phases  respectively  are  decreasing  in  time  due  to  the  wall  stress. 
Figs.  16(c)  and  (d)  show  the  corresponding  radial  velocities.  They  indicate 
the  gas  phase  has  negative  radial  velocity  and  the  liquid  phase  has  positive 
radial  velocity.  That  is,  gas  is  moving  inwardly  and  liquid  is  moving 
outwardly.  Also,  except  near  t = 0,  the  strength  of  these  movements  is 
decreasing.  This  decrease  is  mainly  due  to  the  decreasing  of  the  tangential 
velocity.  Near  t = 0,  the  radial  velocity  is  accelerating  since  the 
analysis  starts  with  no  radial  velocity. 

As  the  result  of  these  radial  movements,  the  volume  fraction 

distribution  is  also  changed  with  time.  Fig.  17  shows  the  gas  volume 

fraction  distributions  of  test  #148  as  functions  of  time.  As  expected,  the 

gas  volume  fraction  is  increasing  at  the  inner  region  and  decreasing  at  the 

outer  region.  However,  because  of  the  decaying  vortex  as  the  result  of  the 

wall  stresses,  the  minimum  gas  volume  fraction  can  only  reach  about  0.18 

even  near  the  wall.  That  is,  for  this  case  a complete  separation  between 

gas  and  liquid  phases  seems  to  be  impossible.  Fig.  18  shows  the  gas  volume 

fraction  distributions  of  test  #150  as  functions  of  time.  Since  a free 

surface,  (I  =0)  was  assumed  for  this  case,  the  vortex  strength  remains 
w 

nearly  constant.  As  the  result  of  the  constant  rotational  strength,  the 
radial  movements  continue  and  the  gas  volume  fraction  approaches  zero  near 
the  wall.  That  is,  there  is  a liquid  layer  formed  near  the  wall  region  for 
this  case. 

Fig.  19(a)  shows  the  liquid  tangential  velocity  of  test  #152.  Near  the 
region  of  injection  (0.85  < r < 0.95),  the  tangential  velocity  is  increasing 
from  zero  very  quickly  at  first  and  then  approaching  an  asymptotic  value. 
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Meanwhile  the  angular  momentum  (or  the  tangential  velocity)  is  transferred 
and  diffused  into  the  inner  region.  Fig.  19(b)  shows  the  liquid  radial 
velocity  of  the  same  test.  Although  the  maximum  tangential  velocity  always 
appears  near  r = 0.85,  the  location  of  the  maximum  radial  velocity  continues 
to  move  inwardly.  This  can  be  explained  more  easily  from  the  time  evolution 
of  the  volume  fraction  distributions  as  shown  in  Fig.  20.  Except  t = 0, 
there  is  an  annular  layer  formed  near  the  wall,  in  which  the  gas  volume 
fraction  is  zero.  This  annular  liquid  layer  grows  with  time,  and  the 
interfacial  radius  of  this  layer  is  moving  inwardly.  Since  a liquid  layer 
contains  little  or  no  gas,  the  associated  radial  movement  will  be  very  small 
even  if  the  flow  has  a large  angular  rotation.  Therefore,  the  location  of 
the  maximum  radial  velocity  is  moving  inwardly  despite  the  fact  that  the 
location  of  the  maximum  tangential  velocity  remains  nearly  unchanged  at  r = 
0.85. 

For  comparison,  the  liquid  velocity  distributions  of  test  #151  for  t = 

4 and  5 are  also  shown  on  Figs.  19(a)  and  (b).  The  conditions  for  #151  are 
same  as  those  of  #152  except  the  injection  rate  is  shut  off  after  t = 3. 

Thus  for  t < 3,  all  the  data  are  exactly  the  same  as  those  of  #152.  Due  to 
the  effect  of  the  viscosity,  the  peak  momentum  will  continue  to  diffuse  and 
transfer  to  the  region  of  lower  momentum,  particularly  at  the  wall  where  a 
momentum  sink  is  located.  Figs.  19(a)  and  (b)  show  that  both  the  tangential 
and  radial  components  of  the  velocity  are  decreasing  with  time  as  the 
injection  is  turned  off  at  t = 3*  However,  the  decrease  of  the  tangential 
velocities  is  much  larger  than  that  of  the  radial  component.  Fig.  20  also 
shows  the  gas  volume  fraction  distributions  of  #151.  Even  the  velocity 
fields  are  quite  different  between  cases  #151  and  #152  especially  for  the 
tangential  component.  The  volume  fraction  distributions  between  these  two 
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cases  are  practically  the  same.  This  indicates  the  final  part  of  the 
injection  (i.e.,  t > 3)  is  very  ineffective. 

Fig.  21  shows  the  gas  volume  fraction  distributions  of  test  #153.  This 
shows  that  a larger  gas  eddy  viscosity  can  produce  a better  separation, 
especially  for  the  region  far  away  from  the  region  of  the  injection.  Fig. 

22  shows  the  gas  volume  fraction  distributions  of  test  #154.  As  expected, 
lower  injection  speed  will  produce  smaller  phase  separation.  However,  the 
effect  is  not  linear.  From  Figs.  (20)  and  (22)  the  time  required  to  produce 
a liquid  layer  of  15%  of  radius,  is  about  1 for  #152  and  about  3 for  #154. 
That  is,  a factor  of  2 reduction  on  the  injection  speed  will  require  a 
factor  of  3 increase  on  time. 

To  study  fluid  dynamical  behavior  of  these  two-phase  fluids,  it  is 
useful  to  know  the  velocity  fields.  To  design  a practical  tank  system,  it 
is  important  to  know  the  pressure  distribution  in  the  tank.  Figs.  (23)  and 
(24)  show  the  pressure  distributions  of  tests  #148  and  #152,  respectively. 
The  pressure  distributions  of  test  #151  at  t = 4 and  5 are  also  shown  in 
Fig.  24  for  comparison. 

As  expected  the  data  show  the  maximum  pressure  occurs  at  the  tank  wall 
for  all  cases.  For  the  case  of  initially  pure  rotation  (#148),  the  pressure 
field  is  decreasing  as  the  rotational  speed  reduces.  The  slightly  different 
character  that  appears  in  the  pressure  distribution  at  t = 0 is  due  to  the 
assumption  of  no  initial  radial  velocity.  The  initial  pressure  field 
essentially  is  that  of  the  centrifugal  forces  due  to  the  pure  rotation.  As 
time  goes  by,  the  radial  velocity  increases  and  the  pressure  field  is  thus 
adjusted  accordingly.  In  the  injection  cases,  the  pressure  field  gradually 
increases  as  the  injection  continues.  As  soon  as  the  injection  is  stopped, 
the  pressure  field  is  changes  as  shown  by  #151  in  Fig.  24. 
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As  discussed  earlier  regarding  the  operation  concept,  the  time 
evolution  of  vortex  strength  and  the  degree  of  phase  separation  are  the 
important  operation  parameters.  The  definitions  of  these  parameters  are 
arbitrary  and  depend  on  various  situations  or  requirements.  Some  of  these 
quantities  are  defined  below: 


a)  Vortex  Strength  r,  the  total  angular  momentum  per  unit  volume: 

f 2 

2 2 a.V..r  dr 

J i 9i 


1 VA/  = 1 P,  ? ? 

i-1  9 i-1  (R2-  R, 2 ) 


b)  Separation  Strength  Q , the  averaged  radial  flow  rate  per  unit  volume: 


Q„  = 


j 2lrra2Vr2dr 


2 2 

ir(R  -Rp 


dr 


2 J a2Vr2rdr 

(R2-R1)(  R-R.^ 


This  indicates  the  rapidity  with  which  the  liquid  phase  moves  to  the 
wall  region  and  the  gas  phase  moves  to  the  inner  region  of  the  circular 


tank.  Thus  a larger  means  a fast  phase  separation. 


c)  Liquid  Purity,  (r^rp,  the  averaged  liquid  volume  fraction  between 


the  radii  of  r^  and  r^: 


ci  2 ( ^ ^ ~ 


fr2 

2 Jr^o^rdr 

~2  2~ 

r2~  V 


For  the  data  presented  below,  r^  = 0.8  and  r^  = 1 are  chosen.  This  is 
the  region  where  the  withdrawal  port  will  most  likely  be  located. 
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d)  Separation  Index,  SPI. 


SPI  = 


VR1 

^2s~^ls 


where 


R - i-S 

v 


a,  rdA 


jakdA 


is  the  center  of  gravity  of  phase  k,  and  R is  the  center  of  gravity  of 

f\  O 

phase  k when  the  two  phases  are  completely  separated. 

When  two  phases  are  completely  separated,  the  liquid  forms  an  annulus 
against  the  wall  of  the  tank  around  a core  of  gas  with  a contact  radius  of 

2 f 1/2 

R = (R  - 2 a_rdr)  . 
c J 2 


Thus 


R.  = 2/3  (R  3-R  3)/(R  2-R  2) 
Is  c 1 c 1 


R0  = 2/3  (R3-R  3)/(R2-R  2) 
2s  c c 


Thus,  SIP  = 0 indicates  that  the  two  phases  are  uniformly  mixed  together, 
SPI  = 1 means  two  phases  are  completely  separated,  and  SPI  < 0 means  more 
gas  is  in  the  outer  region  than  in  the  inner  region, 
e)  Mixing  Index,  MXI 

a A ! ot  a-  dA 

MXI  = = j ‘ — - — t 

[ a.dA][  a_dA] 

a10t2  J 1 J 2 


This  is  another  quantity  to  identify  the  degree  of  separation.  When  mixing 
index  equal  to  1 (MXI  = 1)  it  means  two  phases  are  completely  mixed  and  the 
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volume  fractions  are  uniformly  distributed.  When  mixing  index  decreases  to 
0 (MXI  = 0)  it  means  there  is  no  mixing  between  the  two  or  they  are 
separated  completely. 

Thus,  some  of  the  fluid  dynamical  behavior  of  the  two  phase  fluids  in 
the  tank  can  be  identified.  These  are  the  vortex  strength  r,  the  separation 

strength  Q , the  liquid  purity  a^CO.S,  1),  the  separation  index  SPI,  and  the 
mixing  index  MXI.  These  are  plotted  versus  time  for  various  tests  as  shown 
on  Figs.  (25)— (29) . In  order  to  present  these  using  the  same  ordinate 

scale,  the  quantity  50Q^  (rather  than  Q ) was  plotted.  For  test  # 148,  Fig. 
25  shows  that  the  vortex  strength  is  decreasing,  and  the  separation  strength 
reaches  a maximum  at  t * 0.5  and  then  decreases  monotonically . It  also 

shows  both  liquid  purity  (0.8,  1)  and  separation  index  SPI  are  continuing 
to  improve  even  though  the  vortex  strength  is  decreasing.  However,  because 
the  separation  strength  is  small  and  decreasing,  the  mixing  index  maintains 
a very  high  value  ( ~ 1)  and  the  two  phases  are  still  mixed  very  well  even  at 
t = 5.  (MXI  = 0.95  and  SPI  = 0.27  at  t - 5). 

If  the  effect  of  the  wall  stress  is  removed  as  shown  in  test  #150,  the 
situation  is  somewhat  different.  Fig.  26  indicates  the  vortex  strength 

remains  nearly  unchanged.  Liquid  purity  ap  (0.8,  1),  separation  strength 
and  separation  index,  SPI  are  all  improved  much  more  as  compared  to  Fig.  25. 
Although  the  two  phases  are  not  completely  separated  as  indicated  by  both 

the  values  of  SPI  and  MXI,  the  liquid  purity  ap  has  nearly  reached  the  value 
of  1 when  t > 4.  This  suggests  that,  even  when  the  two  phases  are  not 
completely  separated,  the  liquid  transfer  could  be  started  provided  the 
withdrawal  port  is  located  near  the  wall.  Fig.  27  shows  the  same  quantities 
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for  tests  #152.  For  t > 3,  the  values  for  test  #151  are  also  shown  on  the 


same  figure.  The  figure  shows  there  is  little  or  no  difference  on  the 
degree  of  separation  between  tests  #151  and  #152,  although  the  corresponding 
vortex  strengths  are  quite  different  for  t > 3.  This  indicates  again  that 
the  continuation  of  the  injection  of  #152  for  t > 3 is  not  very  effective. 

Fig.  28  shows  the  results  of  test  #153.  As  compared  to  Fig.  27,  the 
larger  eddy  viscosity  factor  improves  the  degree  of  separation.  This  is 
partly  due  to  the  more  effective  injection  as  the  result  of  the  larger 
diffusion  of  the  injected  angular  momentum  produced  by  the  larger  effective 

viscosity.  Despite  a better  degree  of  separation,  the  values  of  (0.8,  1) 
are  practically  the  same.  That  is,  the  eddy  viscosity  does  not  noticeably 
affect  the  flow  dynamics  in  the  region  of  0.8  < r < 1 near  which  the 

injection  was  made.  The  value  of  a 2 (0.8,  1)  indicates  the  liquid  layer  at 
the  tank  wall  reaches  20$  of  radii  thick  at  t = 1 for  both  cases. 

Fig.  29  shows  the  degree  of  phase  separation  for  the  lower  injection 
case,  test  #154.  As  expected,  the  degree  of  phase  separation  increases 

slower  than  in  the  higher  injection  cases.  Nevertheless,  liquid  purity  otp 
(0.8,  1)  reaches  1 at  t 3 as  opposed  to  at  t * 1 for  the  higher  injection 
cases . 

3 . 7 Conclusions 

A vortex-induced,  liquid  handling  process  adapted  to  two-phase  fluids 
in  zero  gravity  has  been  analyzed  for  a selected  range  of  fluid  and  flow 
parameters.  This  transfer  process  has  been  divided  in  two  stages:  the 
initial  "spin-up"  stage  and  the  liquid  pump-out  stage.  The  initial  spin-up 
stage  is  modeled  in  the  circular  tank  via  tangential  injection  techniques. 
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The  model  is  based  on  a two-phase,  two-fluid  continuum  and  it  incorporates 
several  interactions  between  phases  - namely  fluid  drag  and  pertinent 
virtual  mass.  The  model  provides  a mathematical  framework  containing  the 
essential  physical  mechanisms  on  vortex-induced  phase  separation. 

A computer  program  using  a FORTRAN  77  based  computer  code  and  HP-1000 
minicomputer  was  developed.  The  program  is  flexible  and  accepts  various 
input  parameters  for  different  flow  conditions.  It  can  also  be  modified 
easily  to  include  additional  interaction  effects.  This  program  can  be  used 
to  study  the  fluid  dynamical  behavior  of  two-phase  fluids  in  the  sender 
tank.  It  provides  a quick/easy  sensitivity  test  on  various  parameters  and 
thus  provides  the  guidance  for  the  design  and  use  of  actual  physical  systems 
for  handling  two-phase  fluids  in  space. 

Several  examples  are  given  to  demonstrate  the  products  of  the  program. 
The  interesting  and  plausible  results  seem  to  indicate  the  simplified  two- 
phase  two-fluid  model  is  a useful  one.  However,  detailed  evaluation  of  the 
result  is  only  possible  with  the  aid  of  accurate  local  measurements.  It  is 
therefore  concluded  that  this  model  should  receive  further  attention  in  the 
form  of  experimental  verifications.  These  verifications  should  be  devised 
to  check  both  the  formulations  selected  for  phase  interactions  as  well  as 
the  results  predicted. 

Acknowl edgments 

We  would  gratefully  acknowledge  the  support  received  from  NASA’s 
Kennedy  Space  Center  to  carry  out  the  model-development  work  described 
herein.  Specific  thanks  go  to  Mr.  Frank  Howard  whose  involvements  and 
inputs  on  this  project  have  been  most  helpful.  We  would  also  like  to  thank 
Drs.  James  M.  McMichael  and  L.  Patrick  Purtell  for  their  invaluable 


53 


contributions  on  the  initiation  and  the  early  stage  of  this  project.  In 
addition,  thanks  are  extended  to  Dr.  George  E.  Mattingly  for  his 
encouragement  and  technical  discussions;  to  Mr.  Mike  Hiza  for  his 
involvement  on  the  project  initiation;  and  to  Ms.  Susan  Johnson  for  typing 
the  manuscript. 


54 


References 


1.  Stark,  J.  A.  (1976)  "Low-G  Fluid  Transfer  Technology  Study",  NASA  CR- 
135020. 

2.  Wallis,  G.  B.  (1982)  "Review  - Theoretical  Models  of  Gas-Liquid  Flows", 
ASME  J.  Fluid  Eng.  104,  279. 

3.  Biesheuvel,  A.  and  Van  Wijngaarden,  L.  (1984)  "Two-Phase  Flow  Equations 
for  a Dilute  Dispersion  of  Gas  Bubbles  in  Liquid",  J.  Fluid  Mech.  Vol . 
148,  301-318. 

4.  Delhaye , J.  M.  (1981)  "Basic  Equations  for  Two-Phase  Flow  Modeling", 
Two-Phase  Flow  and  Heat  Transfer  in  the  Power  and  Process  Industries, 
Bergles,  A.  E.  et  al . , Eds.,  Hemisphere  Publishing  Corporation  & 
McGraw-Hill  Book  Company,  Chapter  2,  40-97. 

5.  Dobran,  F.  (1985)  "Theory  of  Multiphase  Mixtures",  Int . J.  Multiphase 
Flow,  Vol.  11,  no.  1,  1-30. 

6.  Drew,  D.  A.  (1983)  "Mathematical  Modeling  of  Two-Phase  Flow",  Annual 
Review  of  Fluid  Mechanics,  Vol.  15. 

7.  Van  Deemter,  J.  J.  and  Van  Der  Laan , E.  T.  (I960)  "Momentum  and  Energy 
Balances  for  Dispersed  Two-Phase  Flow",  Appl.  Sci.  Res.,  Sec.  A,  Vol. 
10,  102-108. 

8.  Yadigaroglu,  G.  and  Lahey,  R.  T.  (1976)  "On  the  Various  Forms  of  the 
Conservation  Equations  in  Two-Phase  Flow",  Int.  J.  Multiphase  Flow, 

Vol.  2,  477-494. 

9.  Zuber , N.  (1964)  "On  the  Dispersed  Two-Phase  Flow  in  the  Laminar  Flow 
Regime",  Chem.  Eng.  Sci.,  19,  897. 

10.  Van  Der  Welle,  R.  (1981)  "Turbulence  Viscosity  in  Vertical  Adiabatic 
Gas-Liquid  Flow",  Int.  J.  Multiphase  Flow,  Vol.  7,  pp.  461-473. 


55 


11.  Brenner,  H.  (1961)  "The  Oseen  Resistance  of  a Particle  of  Arbitrary- 
Shape",  J.  Fluid  Mech.  11,  604. 

12.  Cook,  T.  L.  and  Harlow,  F.  H.  (1984)  "Virtual  Mass  in  Multiphase  Flow" 
Int.  J.  Multiphase  Flow,  Vol.  10,  No.  6,  691. 

13.  Frankel , N.  A.  and  Acrivos,  A.  (1967)  "On  the  Viscosity  of  a 
Concentrated  Suspension  of  Solid  Spheres",  Chem.  Eng.  Sci.,  22,  847. 

14.  Ishii,  M.  and  Zuber , N.  (1979)  "Drag  Coefficient  and  Relative  Velocity 
in  Bubbly,  Droplet  or  Particulate  Flows",  AIChE  J.  25,  pp.  843-55. 

15.  Lee,  N.  and  Dukler,  A.  E.  (1981)  "A  Stochastic  Model  for  Turbulent 
Diffusion  of  Particles  or  Drops"  AIChE  J.  Vol.  27,  No.  4,  552. 

16.  Magnus,  R.  J.  (1966)  "Analytic  Prediction  of  Pressure  Pulses  Due  to  the 
Rocket  Engine  Starting",  General  Dynamics,  GDC-ERR-AN-984 . 

17.  MacCormack,  R.  W.  (1969)  "The  Effect  of  Viscosity  in  Hypervelocity 
Impact  Cratering",  AIAA  Paper  6 9 —3 54. 

18.  Chima,  R.  V.  and  Johnson,  G.  M.  (1985)  "Efficient  Solution  of  the  Euler 
and  Navier-Stokes  Equations  with  a Vectorized  Multiple-Grid  Algorithm", 
AIAA  Vol.  23,  No.  1. 

19.  Davis,  R.  W.  (1983)  "Finite  Difference  Methods  for  Fluid  Flow", 
Proceedings  of  the  1983  International  Conference  on  Computational 
Techniques  and  Applications",  edited  by  John  Noye  and  Clive  Fletcher, 
pp.  51-69. 

20.  Peyret , R.  and  Taylor,  T.  D.  (1983)  "Computational  Methods  for  Fluid 
Flow"  Spr inger-Verlag , New  York. 


56 


Vortex 

Rotation  Gas  Injection  Gas 


cc 

v 

\ 

\ 

\ LlI 

• « , 

o / 

o 

3? 

O 

3 

O’ 

ii/  Ld 

/ 

rr 

fj  ££ 

i . t ^ 


e 

u 

co 

cn 

u 

oj 

co 

c 

03 

U 

H 

H 

co 

i— i 

a 


<u 

u 

3 

ao 

•H 


57 


58 


Figure  2.  Tangential  Vortex  Model 


Figure  3.  Start-up  Configuration 
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Figure  4.  Pump-  out  Configuration 
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Figure  5.  The  Two-Stage  Spin- up  and  Transfer  Concept 
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Figure  8.  Modeling  of  Drag  Coefficient 
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Figure  10.  Two  Step  Difference  Scheme  (Backward  Predictor  - 
Forward  Corrector  Version) 
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Figure  11.  Velocity  Vector  Distributions  of  Test  #148 
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Figure  13.  Velocity  Vector  Distributions  of  Test  Run  #151 
The  annular  region  between  two  dashed  circles  is 
the  region  of  injection. 
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the  region  of  injection. 
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Figure  15.  Velocity  Vector  Distributions  of  Test  Run  //153 
The  annular  region  between  two  dashed  circles  is 
the  region  of  injection. 
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Figure  16.  Velocity  Distributions  of  Test  #148 
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Figure  18.  Gas  Volume  Fraction  Distributions  of  Test  //150 
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Liquid  Tangential  Velocity  Distributions  of  Test  #152 
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Figure  19(b).  Liquid  Radial  Velocity  Distributions  of  Test  #152 
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Figure  20.  Gas  Volume  Fraction  Ois t rib ut ions  of  Test  #152 
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Figure  21.  Gas  Volume  Fraction  Distributions  of  Test  7/153 
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t’  igure  22*  Gas  Volume  Fraction  Distributions  of  Test  //154 
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Figure  23.  Pressure  Distributions  of  Test  #148 
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Figure  24 . Pressure  Distributions  of  Test  #153 
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Figure  25.  The  Vortex  Strength  and  Quantities  on  the  Degree  of  Phase  Separation  of  Test  // 148 


F i gure  26.  The  Vortex.  Strength  and  Quantities  on  the  Degree  of  Phase  Separation  of  Test  # JL  i>0 
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Figure  27.  The  Vortex  Strength  and  Quantities  on  the  Degree  of  Phase  Separation  of  Test  #152 
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Figure  28.  The  Vortex  Strength  and  Quantities  on  the  Degree  of  Phase  Separation  of  Test  //153 
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Figure  29.  The  Vortex  Strength  and  Quantities  on  the  Degree  of  Phase  Separation  of  Test  //154 
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with  several  phase  interactions  - namely  fluid  drag  and  pertinent  virtual  mass. 

A computer  program  was  developed  to  study  the  fluid  dynamical  behavior  of 
two-phase  fluids  in  the  sender  tank.  Several  examples  are  given  to  demonstrate  the 
products  of  the  program.  The  interested  and  plausible  results  indicate  the 
simplified  two-phase  two-fluid  model  is  a useful  model.  However,  a detailed 
evaluation  of  the  result  is  only  possible  with  the  aid  of  accurate  local  measure- 
ments. These  verifications  should  be  devised  to  check  both  the  formulations 
selected  for  phase  interactions  as  well  as  the  results  predicted. 
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